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Abstract 

The gradient of potential vorticity (PV) is an important quantity because of the way 
PV (denoted as q) tends to accumulate locally in the oceans and atmospheres. Re- 
cent analysis by the authors has shown that the vector quantity B = Vg x V0 for the 
three-dimensional incompressible rotating Euler equations evolves according to the same 
stretching equation as for u the vorticity and B, the magnetic field in magnetohydro- 
dynamics (MHD). The B-vector therefore acts like the vorticity u in Euler's equations 
and the B-ileld in MHD. For example, it allows various analogies, such as stretching dy- 
namics, helicity, superhelicity and cross helicity. In addition, using quaternionic analysis, 
the dynamics of the B-vector naturally allow the construction of an orthonormal frame 
attached to fluid particles; this is designated as a quaternion frame. The alignment 
dynamics of this frame are particularly relevant to the three-axis rotations that particles 
undergo as they traverse regions of a flow when the PV gradient V(7 is large. 
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1 Introduction 

The ideas in this paper weave together two strands of research on the Euler fluid equations 
recently pursued by the authors. Both of these required the use of Ertel's Theorem [Ij - 
always a favourite topic with Raymond - and were discussed at length with him during their 
development. 

The first and latest strand of research concerns the evolution of the gradient of potential 
vorticity (Vg) and the gradient of potential temperature (V^^) [5]. In the case of the Euler 
equations, while both q and are material constants, the evolution of their gradients involves 
the strain and rotation rates of the flow. Physically, understanding the behaviour of Vg in the 
atmosphere and the oceans is of paramount importance because potential vorticity tends to 
accumulate into localised spatial regions (patches) with sharp edges, where the magnitude | Vg| 
is much larger than its average value [3J. In this regard, the divergenceless vector combination 
B = Vg X V^^ is a natural choice, because it leads to an evolution equation identical to that 
for either the vorticity in the incompressible three-dimensional Euler equations, or the magnetic 
field in an electrically conducting fluid. As a consequence, 3 may undergo the same violent 
stretching and twisting associated with the vorticity field in three-dimensional turbulence, or 
with magnetic field lines in magnetohydrodynamics (MHD), particularly if B were to align 
with an eigenvector of the 3x3 strain-rate matrix S associated with the fluid motion. 

The second strand of research involves the use of quaternions in identifying an ortho-normal 
frame attached to fluid particles in an Euler flow and whose dynamics represent the tumbling 
of the particle as it undergoes three-dimensional rotations during its flight [4l[5]. The vector 
B turns out to be an ideal candidate for the construction of this ortho-normal frame. This 
would be applicable to the dynamics of particles travelling through regions of the oceans or 
atmospheres which have a high values of |Vg|. 

1.1 Potential vorticity gradient for the incompressible Euler 

Consider the dimensionless form of the Euler equations for incompressible, stratified and ro- 
tating flow 

Du D d 

— + 2nxu + aok9 = -Vp, _ = — + n-V, (1.1) 

Dt Dt ot 

in which the potential temperature 9{x, t) evolves passively according to 

DO , , 

In (II. ip . the vectors fi and k are the rotation rate and vertical direction, respectively, and the 
scalar is a dimensionless constant. Information about V^^ would be needed to discuss how 
0{x, t) might accumulate into large local concentrations. This is best studied in the context 
of potential vorticity defined by 

q := u^rot ■ Ve (1.3) 

where ojrot is defined as ujrot := ^ + 217 and uj := curlw denotes the fluid vorticity. Ertel's 
theorem says that the material time derivative D/Dt and the vector field iVrot ■ V operating 
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on a scalar function commute with each other under an Euler flow |T] 

- LJrot ■ ] -^6 + UJrot " ^ TiT • (1-4) 



Dt V Dt J \Dt J 

When UJrot obeys the incompressible Euler equations (V"*" = {dy, —dx, 0)) 

D UJrot 



Dt 

then g is a also material constant because 

Dq 



UJrot ■ Vu - aoV^e (1.5) 




^ , • 1- X V - const 

: const 



Figure 1 : The vector B = Vf/ x is tangent to the curve defined by the intersection of the two 
surfaces q = const and 9 = const in three dimensions. 

To achieve this, the following divergence-free flux vector is constructed [2] 

B = VQ{q) xVe, (1.7) 

where Q{q) is any smooth function of g. The vector B, as in Figure 1, lies along the intersection 
of iso-surfaces of q and 9 and could be thought of as pointing along the tangent to an iso-PV 
curve on a level set of temperature, or vice versa. Thus, the PV function Q{q) is the stream 
function for the flux vector (B) on a level set of potential temperature (9). This observation 
results in the remarkably simple evolution equation 

dB 

— - curl (n X B) = -V(gQ'divn) x . (1.8) 

The cross-product combination is special : an appendix in [5] contains a proof of (11.81) per- 
formed using differential geometry and also by conventional vector analysis. It follows from 
(II. 8p that B satisfies the stretching relation 

DB 

— - = B ■ Vu - B div u - V(gg'div u)xV9, (1.9) 
whose properties will be discussed in the next sectior0. 



^In |2] the original references have been discussed ([6] |7] [8]) in which (jl.Sp had been derived for the 
incompressible Euler equations where the right hand side turned out to be zero because divit = 0. In [2\ the 
argument is extended to both the Navier-Stokes and hydrostatic viscous primitive equations. 
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1.2 Stretching, helicity, superhelicity and cross helicity 

In the incompressible case where divn = 0, equation (11.90 simplifies so that the divergenceless 
vector 3 satisfies the same stretching equation as that for the vorticity u ; namely 

DB 

-— = BVu with divB = 0. (1.10) 

Dt ^ ' 

This leads immediately to 

^|B|2 = e.^B = A(^)|B|2, (1.11) 

where \^^\x, t) is an estimate for an eigenvalue of the rate of strain matrix S and lies within 
its spectrum. Alignment of S with a positive (negative) eigenvector of S may produce violent 
growth (decay) thus reproducing the stretching mechanism that produce the very large vorticity 
intensities that can develop locally in the early and intermediate stages of turbulence. 

Moffatt's analogy between vorticity and magnetic field [9], and his detailed discussion of 
the topology of magnetic field lines, is based on the concept of helicity which for us requires 
the existence of a vector potential A. defined by 

A=\{QVe-dVQ)+Vi), (1.12) 

in terms of the dynamical quantities Q, 9 and a potential ip. Then the helicity H, defined by 



H= A-BdV=(b ^pB-hdS, (1.13) 

Jv JdV 

measures the number of linkages of the field lines of B with themselves. The time derivative 
of helicity under the flow of the Euler equations is given by 



dH 
~dt 



I -{A-B)u-h+ (^)B-h 

JdV L \ Dt / 



dS, (1.14) 



which would vanish for either homogeneous or periodic boundary conditions. For the Euler 
equations, n ■ n = is imposed on a fixed boundary. However, H ^ and dH/dt ^ would 
be possible for topography in which B • n 7^ 0, so the boundaries play the only role in allowing 
linkages in the B-field as there is no other source of helicity. 

Hide's intriguing concept of a super-helicity, which measures the linkages of the field lines 
of JT" := curIB with itself, may be introduced for the jT-vector as in for MHD [TU]. Super- 
helicity is defined as 

S = [ B-JdV. (1.15) 
Jv 

After a short computation, the super-helicity dynamics for the ^/-vector comes out to be 

dS 



2B-cur\\uxB)dV + (t {u ■ J) B ■ h - {J ■ B) u ■ h dS , (1.16) 

V JdV 



which, unlike the helicity H, has both volume and surface sources. Likewise the cross helicity 
for the B-vector can be introduced 

C = [ u-BdV (1.17) 
Jv 
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in analogy with the corresponding quantity in MHD. Another short computation produces the 
dynamics of the cross helicity (R is the vector potential such that curl R = 20) 

^ = -ao J {eB-k)dV + ^ (^-p + u-R+ ^u'^^B-hdS , (1.18) 

which again has both volume and surface sources. 

1.3 Quaternions and an attached orthonormal frame 




trajectory of tracer particle 



Figure 2 : The solid curves represent characteristic curves B = dx/ds (s is arc length) to which 
B is a unit tangent vector. The quaternion-frame orientation {B, x^, B x x^) is shown at the two 
space-time points {xi,ti) to {x2,t2); note that this is not the Frenet-frame corresponding to the 
particle path but to curves B = dx/ds. The dotted line represents the tracer particle (•) path. 

More than one hundred and fifty years ago William Rowan Hamilton invented quaternions 
as a means of representing a composition of rotations. For most of this period they have 
been under-appreciated, yet they have recently undergone a spectacular renaissance due to 
their efficacy in certain applications in avoiding the difficulties incurred at the north and south 
pole^ when Euler angles are used in computing the dynamics of objects undergoing three-axis 
rotations [llj. In particular, quaternions now lie at the heart of many modern inertial guidance 
systems where tracking the paths and the orientation of aircraft and satellites is of importance 
[12] . They are also used in the graphics community to control the orientation of tumbling 
objects in computer animations 

A natural question is whether quaternions are useful in tracking the angular velocity and 
orientation of Lagrangian particles in fluid dynamics. Experiments in turbulent flows have now 
reached the stage where the trajectories of tracer particles can be detected at high Reynolds 
numbers [131 [IS [l5l [T6J. Numerical differentiation of these trajectories gives information 
about the Lagrangian velocity and acceleration of the particles and also the curvature of the 
particle paths [16J. 

Conventional practice has been to consider the Frenet-frame of a trajectory. This consists 
of a unit tangent vector, a normal and a bi-normal, which are used to represent the pitch, yaw 
and roll of the motion. While the Frenet-frame describes the path, it ignores the rotational 
dynamics of the particle. To account for this, another ortho-normal frame associated with the 

^Computations with Euler angles often suffer from "gyro-lock" because of singularities at the poles of the 
spherical angular coordinate system where the azimuthal angle is undefined. 
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motion of a Lagrangian fluid particle - designated the quaternion-frame - has been introduced 
by the authors [4]. This frame moves with the particles, but its evolution derives from the 
fluid equations of motion. In the context of the incompressible Euler equations, the idea of 
the quaternion frame depends on the existence of vectors for which there exists a Lagrangian 
equation of motiorlfj. The natural triplet is {uj, X, ci; x where x = ^ ^ Su: - see [U [5]. 
For the incompressible Euler equations the natural candidate for the B-stretching equation 
( inU|l is the triplet 



{B, Xa, ^ X Xa} where Xa = ^ x « (1-19) 

where the a-label has its origins in the definition a := B ■ Vu. The Lagrangian equations 
of motion for this triplet is derivable through a quaternionic formulation of the Lagrangian 
equation of motion for B given in (ll.lOp . 



2 Quaternions, rigid body rotations and their properties 

The material in this section provides the reader with a definition of quaternions, together with 
a precis of their multiplication rules and properties. 

The literature on rotations in rigid body mechanics is replete with explicit formulae relating 
the Euler angles and what are called the Cayley-Klein parameters of a rotation [IB]. The 
complicated inter-relations that are unavoidable when Euler angle formulae are used can be 
avoided when quaternions are used [13 [20] : for a more modern context see Holm [211 122] 
and Marsden and Ratiu [23j . 

In terms of any scalar p and any 3-vector q, the 4-vector quaternion q = [p, q] is defined 
as (Gothic fonts denote quaternions) 

3 

^ = [p, q]=p\-'^(ii(^i, (2.1) 

1=1 

where {ai, (72, 0-3} are the three Pauli spin-matrices defined by 

-^C-O' 

and I is the 2x2 unit matrix. A rule (denoted as ®) for multiplication of quaternions can be 
derived from the multiplication rule for the Pauli matrices a^aj = —5ij\ — eijkO'k, as 

qi ® q2 = [piP2 - Qi ■ 92, Piq2 + P2qi + qi X 92] ■ (2.3) 

Thus, the multiplication of quaternions is associative, but it is not commutative. 

Let p = [p, q] be a unit quaternion, satisfying p^ + = 1. Its inverse p* = [p, —q] 
satisfies p ® p* = [p^ + g^, 0] = [1,0]. A pure quaternion has a zero scalar entry, such as 

■^The eigenvectors of the rate of strain matrix S are excluded for this reason. 
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X = [0, r]. Hamilton called his pure quaternions vectors. Pure quaternions transform among 
themselves as r = [0, r] — )■ 9^ = [0, R] under 

fH = p®r®p*. (2.4) 

This associative product can be written as 

= p ® r ® p* = [0, (/ - g^)r + 2p{q x r) + 2q{r ■ q)] . (2.5) 
Choosing p = ± cos ^9 and q = ±nsm \6, where n is the unit normal to r, we find that 

01 = p @r ®p* = [0, rcos0 + (n X r)sin^] = 0(0, n)r, (2.6) 

where 

p = ±[cosi6', nsini6'] . (2.7) 

Equation (12.61) is the Euler-Rodrigues formula for the rotation 0(9, h) by an angle 9 of the 
3-vector r about its normal n and the quantities 9 ,fi are the Euler parameters. The elements 
of the unit quaternion p are the Cayley-Klein parameters which are related to the Euler angles 
[18J . When p is time-dependent, the Euler-Rodrigues formula in (12.61) is 

^H(t)=p@r®p* r = p* ® ^H(t) ®p. (2.8) 

It is necessary to use the property of the pure quaternion 91* = —D\ to obtain the time 
derivative of 91 

91(t) = (^®p*)®91-((^®p*)®91)*. (2.9) 

The quaternion p = [p, q] is of unit length and so pp + qq = 0, which means that p ® p* is 
also a pure quaternion 

^®p* = [0, |f7o(t)]. (2.10) 
The 3-vector entry in (I2.10p defines the angular frequency ftoit) as CIq = 2{—pq + qp — q x q) 



thereby giving the well-known formula for the rotation of a rigid body 

R=noxR. (2.11) 

For a Lagrangian particle, the equivalent of Qq is the Darboux vector Da in Theorem [T] in the 
next section. 

3 An ortho-normal frame and particle trajectories 

Having set the scene in ^ by describing some of the essential properties of quaternions, it is 
now time to apply them to the Lagrangian relation (11.101) between the two vectors B and a 

^ = a=BVu. (3.1) 

Dt ^ ' 
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It will turn out below that a knowledge of Da/ Dt is needed. Ertel's Theorem is applicable 
and the result becomes a version of Ohkitani's relation [17j 

^(B.V.) ^ S.V(^) (3,2) 

= -PB - B ■ V (2ri X w + a^ke) (3.3) 

where the Hessian matrix of the pressure p is defined as 

P (3 4) 

Thus we can define a vector b such that 

^ = b:= -PB-B -V {2nx u + aok9) . (3.5) 

Through the multiplication rule in (I2.3p quaternions appear in the decomposition of the 3- 

vector a into parts parallel and perpendicular to another vector, which we choose to be B. 
This decomposition is expressed as 

a = + Xa X B = K, Xa] ® [0, B] , (3.6) 

where the scalar and 3-vector Xa ^''^ defined as 

aa = B-\B-a), Xa = B-\Bxa). (3.7) 

Equation (13. 6p thus shows that the quaternionic product is summoned in naturally0 It is now 
easily seen that is the growth rate of the scalar magnitude {B = \B\) which obeys 

^ = aaB, (3.8) 

Dt ^ ^ 

while Xa' the swing rate of the unit tangent vector B = BB^^, satisfies 

DB 

'Dt 

Now define the two quaternions 



X,xB. (3.9) 



qa=[aa,Xa], ^ = [0, B] , (3.10) 

so (13. ip can automatically be re-written in the quaternion form 

D^ 

-^ = qa®^. (3.11) 



''With reference to fj2] the Cayley-Klein parameters of the quaternion q = [a, x] are 



Q,2 ^ -j^2 ' Q,2 _|_ -^2 
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Moreover, because of (13. Sp . exactly as for q^, a quaternion q?, can be defined which is based 
on the variables 

ai, = B-\B-b), Xb = B-\Bxb), (3.12) 

where 

^b = Xbl ■ (3-13) 
The 3-vector b = Da/Dt admits a decomposition similar to that for a as in (13.61) 

^ = [0, b] = [0, + X6 X B] = q, d) !B . (3.14) 
Using the associativity property, compatibility of (I3.14P and (13. lip implies that {B = \B\ ^ 0) 

^ + qa®qa-qb^ ®23 = 0. (3.15) 
This establishes a Riccati relation between q^ and q^ 

Dt 

with components 



+ qa®qa = q6, (3.16) 



-^K, Xa] + ["a - Xl, '^aaXa] = [«& , Xb] , (3.17) 

where Xa = \Xa\- There follows a Theorem on a Lagrangian particle undergoing fluid motion 
that is equivalent to the well-known formula (I2.1ip for a rigid body undergoing rotation about 
its center of mass : 

Theorem 1 (Alignment dynamics) The ortho-normal quaternion-frame (S, x^, B x 
has Lagrangian time derivatives 



Dt 
DjB X ica) 
Dt 

DXa 



DB 

-- T>aXB, (3.18) 



X (B X Xa) , (3.19) 
^aXXa, (3.20) 



Dt 

where the Darboux angular velocity vector "Da is defined as 

'Da = Xa+—B, Cb = B-{XaXXb), (3-21) 

Xa 

and the quantities [aa, Xa] ^^'^ ['^b, Xt] ^''^ defined in ^3. 7\) and ^3. 12\) . 

Remark: The Darboux vector X>a is driven by the 3-vector b = Da/Dt which sits in in 
(I3.2ip . The analogy with rigid body rotation expressed in (12. lip is clear. 
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Proof: To find an expression for the Lagrangian time derivatives of the components of the 
frame (S, x^, B x requires the derivative of Xa- To find this, it is necessary to use the 
fact that the 3-vector h can be expressed in this ortho-normal frame as the linear combination 

B-^b = a,B + c^Xa + d,{18 X Xa) ■ (3-22) 

where Cf, is defined in (I3.2ip and db = — (Xa " Xb)- The 3-vector product Xb = B^^{B x b) 
yields 

Xb = CbiBxXa)-dbXa- (3.23) 
When split into components, equation (13.171) becomes 

Daa 



and 



Dt 

DXa 

Dt 



2 2 I 



-2aaXa + Xb ■ 



(3.24) 
(3.25) 



A little more working gives the alignment dynamics in equations (I3.18p - (l3.2ip . 



4 Conclusions 

When the quaternion approach to rotations outlined in ^is applied to the Euler equations it 
demonstrates that quaternions are a natural way of calculating the orientation of Lagrangian 
particles in motion through the concept of ortho-normal quaternion-frames attached to each 
particle. In this particular context, where the S-field is the vector that helps us understand the 
evolution of Vg and V6', knowledge of the quartet of 3-vectors (w, B, a, b) is sufficient for 
the application of Theorem [H The complexity of the 3D Euler equations comes through the 
ortho-normal dynamics via the pressure field. In the present state of knowledge, the projection 
P3 that is part of b would have to be found by computational means. 

A natural question is whether these ideas can be applied to the Navier-Stokes equations? 
It turns out that the evolution equation for B in the incompressible case is the same as that 
in (11.91) with u replaced by W |2|. W is a new transport velocity calculated using the method 
of Haynes and Mclntyre [24] 

q{TU-u) = -{ [Re-^/\u + 2Vtxu + a^kO] xVe + {aReY^ujM] (4.1) 

The quaternion procedure can only be pursued to a certain point for the Navier-Stokes equa- 
tions, after which a difficulty appears when DU/Dt is needed, and we have no proper knowl- 
edge of this. This is consistent with the objections that W is not a genuine pliysical velocity 
but merely a mathematical construction [251121]. For the Euler equations we know that this is 
the point where the Hessian matrix of the pressure is introduced in equation (13.20 . Moreover, 
the stretching relation (II. 9p in this case becomes 
DB 

-— = B-VU-Bd\yU-V{qQ'd\vU) x V^, (4.2) 

and divW ^ 0, which allows richer and potentially more singular alignment dynamics than 
those for the incompressible case discussed here in equations (I3.18p -( l?.2ip . As well as the 



Navier-Stokes alluded to above, the case of the hydrostatic primitive equations have been 
discussed in this context in [2]. 
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